Что за случайности и как их генерировать

Филипп

2018-05-31

Чем будем заниматься

  • посмотрим на то, что уже знаем под новым углом

  • немного углубим знания

  • будем смотреть презы, писать код, болтать, решать задачки на доске

  • попробуем быть поближе к реальности

В чём будем заниматься

  • конечно же в R

  • R очень лёгок в освоении

  • одна из самыз красивых визуализаций данных

  • очень хорош в работе со статистикой

  • много готовых пакетов и большое комьюнити

  • используется многими компаниями в работе, вот устаревший мини-список

Про активности

  • не делаешь сам \(\Rightarrow\) никогда не научишься \(\Rightarrow\) три большие домашки

  • домашки делаются в командах по три человека, баллы ставятся на команду

  • для зачёта надо набрать …..

  • две команды, набравшие наибольшее количество баллов свожу в офис Яндекса почилить на крыше под винишко

  • в каждой домашке много рыбёшек и большой кит, можно решать только то, к чему лежит душа

Теория вероятностей

  • мы изучали вероятности

  • мы изучали события

  • мы изучали случайные величины и их распределения

  • ЗАЧЕМ МЫ ЭТО ДЕЛАЛИ?

Случайные ли величины?

Демон Лапласа

Байесовский взгляд на вероятность

  • Лаплас: детерминизм, мы могли бы идеально прогнозировать вселенную, если бы измерили точное положение каждого атома. Издержки этого огромны.

  • Между совершенством природы и несовершенством человеческого познания огромный разрыв.

  • Неопределённость — результат этого разрыва. Случайность — это результат нашего невежества, а вероятность способ это невежество измерить.

Частотный взгляд на вероятность

  • Наука не может рассматривать вероятность как субъективную меру невежества

  • Можно оценивать вероятности только тех событий, которые происходят более одного раза

  • Вопрос «Какова вероятность, что кандидат N победит на выборах?» не имеет ответа, так как событие уникально и не обладает «частотой».

Что мы имеем

  • Сундук — различные процесс порождения данных. Теория вероятностей изучает этот сундук. В реальности мы его не видим.

  • В реальном мире мы видим выборки, которые сундук выплёвывает на нас. Математическая статистика изучает испражнения сундука и по ним пытается восстановить его внутренности.

Хлам из сундука

Чем плюётся сундук

  • Выборки, по которым можно попробовать восстановить его внутренности

Базовые теоремы

  • Все манипуляции по восстановлению внутренностей сундука по его испражнениям позволяет делать ряд теорем

  • Некоторые из них вы уже знаете, что это за теоремы?

Закон больших чисел

Закон больших чисел

  • ЗБЧ утверждает, что среднее арифметическое большого числа похожих случайных величин «стабилизируется» с рочтом их числа

  • Как бы сильно случайные величины не отклонялись от своего среднего значения, эти отклонения взаимно гасятся и среднее арифметическое приближается к постоянной величине

\[ \frac{X_1 + \ldots + X_n}{n} - \frac{E(X_1) + \ldots + E(X_n)}{n} \overset{p}{\to} 0 \]

  • ЗбЧ сформулировано довольно много: Чебышёва, Бернулли, Хинчина и тд

Закон больших чисел

Вопрос про больницы

  • Есть две больницы: маленькая и большая
  • В обеих принимают роды, выяснилось, что в одной из больниц оценка вероятности появления мальчика составила \(0.7\)

  • В какой из больниц это произошло и почему?

Вопрос про больницы

  • Есть две больницы: маленькая и большая
  • Скорее всего, это произошло в маленькой больнице. При малых объёмах выборки вероятность отклониться от \(0.5\) больше. Именно это нашёптывает нам ЗБЧ.

Слабая форма ЗБЧ (Пафнутий Чебышёв)

Пусть \(X_1, \ldots, X_n\) попарно независимые и одинаково распределённые случайные величины с конечным вторым моментом, \(E(X_i^2) < \infty\), тогда имеет место сходимость

\[ \frac{X_1 + \ldots + X_n}{n} \overset{p}{\to} E(X_1) \]

Сильная форма ЗБЧ

Пусть \(X_1, \ldots, X_n\) попарно независимые и одинаково распределённые случайные величины с конечным вторым моментом, \(E(X_i^2) < \infty\), причём \(\sum_{k=1}^{\infty} \frac{Var(X_k)}{k^2} < \infty\), тогда имеет место сходимость

\[ \frac{X_1 + \ldots + X_n}{n} \overset{п.н.}{\to} E(X_1) \]

Страховые компании

  • Примерно с \(1600\)-х годов ЗБЧ позволяет зарабатывать деньги на страховании

  • Упражнение: для \(25\)-летней девушки вероятность прожить ещё год составляет \(0.9\). Страховка в год стоить \(1000\) рублей при взносе в \(110\) рублей. Какой будет средняя прибыль компании с одной страховки?

Страховые компании

  • Упражнение: для 25-летней девушки вероятность прожить ещё год составляет \(0.9\). Страховка в год стоить \(1000\) рублей при взносе в \(110\) рублей. Какой будет средняя прибыль компании с одной страховки?

  • \(X_i\) — прибыль с одной страховки

\(X_i\) \(110\) \(-890\)
\(P(\ldots)\) \(0.9\) \(0.1\)

Средняя прибыль компании составит \(\frac{1}{n} \sum X_i\). По закону больших чисел:

\[ \frac{1}{n} \sum_{i=1}^n X_i \overset{p}{\to} E(X_1) = 0.9 \cdot 110 - 0.1 \cdot 890 = 10 \]

Генерации разных штук

  • Другой ништяк, который нам разрешает ЗБЧ — генерация случайных величин для оценки разных математических ожиданий и тп

  • Не можешь посчитать? Сгенерируй!

  • Обычно такие генерации называют методом Монте-Карло

  • Метод Монте-Карло это общее название группы численных методов, основанных на получении большого числа реализаций случайного процесса, который формируется таким образом, чтобы его вероятностные характеристики совпадали с аналогичными величинами решаемой задачи

Генерации в R

Генерации в R

  • Чтобы сварить в R любую случайную величину, нужно знать четыре буквы: r, d, p и q

  • rnorm эта команда сгенерирует выборку из нормального распределения
  • dnorm эта команда вычислит значение плотности в указанной точке
  • pnorm эта команда находит вероятность
  • qnorm эта команда находит квантили

Генерируем нормальную выборку

Хочу сгенерировать нормальную случайную величину

\[ f(x) = \frac{1}{\sigma \sqrt{2 \pi}} \cdot e^{-\frac{(x - \mu)^2}{2 \sigma^2}} \]

x <- rnorm(5, mean=5, sd=3)  # если дисперсия 9, 
x                            # то стандартное отклонение 3
## [1] 5.5863833 2.0927152 5.2230942 0.9752592 5.9821973

Хочу узнать \(f(3)\)

dnorm(3, mean=5, sd=3)                        
## [1] 0.1064827

Генерируем нормальную выборку

Хочу узнать \(F(3) = P(X \le 3)\)

pnorm(3, mean=5, sd=3)                     
## [1] 0.2524925

Генерируем нормальную выборку

Хочу узнать \(F(3)\)

\[ F(3) = P(X < 3) = \int_{-\infty}^3 f(x)dx \]

pnorm(3, mean=5, sd=3)                     
## [1] 0.2524925

Хочу узнать \(P(4 < X < 9)\) - ?

Генерируем нормальную выборку

\[ F(3) = P(X < 3) = \int_{-\infty}^3 f(x)dx \]

pnorm(3, mean=5, sd=3)                     
## [1] 0.2524925

\[ P(4 < X < 9) = \int_4^9 f(x) dx = F(9) - F(4) \]

pnorm(9, mean=5, sd=3) - pnorm(4, mean=5, sd=3)                   
## [1] 0.5393474

Генерируем нормальную выборку

Квантиль уровня \(\gamma\) это такое число \(q\), что:

\[ P(X < q) = \gamma \]

qnorm(0.7, mean=5, sd=3)                   
## [1] 6.573202
x <- c(0.95, 0.975, 0.995)
qnorm(x, mean = 0, sd = 1)                
## [1] 1.644854 1.959964 2.575829

Генерируем нормальную выборку

Что генерирует команда ниже?

rnorm(10)       
##  [1] -1.4795307  0.4837291 -0.7179522  1.0988786  1.0140746  1.0279011
##  [7] -0.8357659 -1.0851605  1.2081977  1.1808669

Генерируем нормальную выборку

Можно посчитать для случайной величины какие-нибудь характеристики

x <- rnorm(1000, mean=5, sd=3)

mean(x) # среднее

var(x) # выборочная дисперсия

sd(x) # выборочное стандартное отклонение

median(x) # выборочная медиана
## [1] 5.018105
## [1] 8.842121
## [1] 2.97357
## [1] 5.097322

Генерируем нормальную выборку

Можно построить гистограму

library("ggplot2")  # пакет для красивых картинок

x <- rnorm(1000, mean=5, sd=3)

qplot(x) # гистограма 
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Генерируем нормальную выборку

Можно построить функцию плотности

library("ggplot2")  # пакет для красивых картинок

x <- seq(-5, 15, by=0.01)
y <- dnorm(x, mean=5, sd=3)

qplot(x, y, geom="line")

Генерируем нормальную выборку

Можно построить функцию распределения

library("ggplot2")  # пакет для красивых картинок

x <- seq(-5, 15, by=0.01)
y <- pnorm(x, mean=5, sd=3)

qplot(x, y, geom="line")

Зачем всё это надо?

  • Чтобы решать реальные проблемы! Например, пусть \(X \sim N(5,3)\). Как найти \(E \left(\frac{1}{X} \right)\)?

  • Спосбо первый:

\[ E \left(\frac{1}{X} \right) = \int_{-\infty}^{+\infty} \frac{1}{x} \cdot \frac{1}{\sigma \sqrt{2 \pi}} \cdot e^{-\frac{(x - 5)^2}{2 \cdot 3^2}} dx \]

Зачем всё это надо?

  • Чтобы решать реальные проблемы! Например, пусть \(X \sim N(5,3)\). Как найти \(E \left(\frac{1}{X} \right)\)?

  • Способ второй:

n_obs <- 10^6 # число наблюдений 

x <- rnorm(n_obs, mean = 5, sd = 3)

mean(1/x)
## [1] 0.2633082

Зачем всё это надо?

  • Чтобы решать реальные проблемы! Например, пусть \(X \sim N(5,3)\). Как найти \(E \left(\frac{1}{X} \right)\)?

  • Способ второй:

n_obs <- 10^6 # число наблюдений 

x <- rnorm(n_obs, mean = 5, sd = 3)

mean(1/x)
## [1] 0.2697192
  • ЗБЧ разрешает нам так делать!

Зачем всё это надо?

  • \(X_1, X_2, X_3 \sim U[0;2]\), независимые

  • Хотим знать \(P(X_1 + X_2 + X_3^2 > 5)\)

Зачем всё это надо?

  • \(X_1, X_2, X_3 \sim U[0;2]\), независимые

  • Хотим знать \(P(X_1 + X_2 + X_3^2 > 5)\)

n_obs <- 10^6
x_1 <- runif(n_obs, min = 0, max = 2)
x_2 <- runif(n_obs, min = 0, max = 2)
x_3 <- runif(n_obs, min = 0, max = 2)

success <- x_1 + x_2 + x_3^2 > 5

success[1:5]

sum(success) / n_obs
## [1] FALSE  TRUE FALSE FALSE  TRUE
## [1] 0.147194

Зачем всё это надо?

  • \(X_1, X_2, X_3 \sim U[0;2]\), независимые

  • Хотим знать \[P(X_1 + X_2 > 0.8 \mid X_3 < 0.1)\]

Зачем всё это надо?

  • \(X_1, X_2, X_3 \sim U[0;2]\), независимые

  • Хотим знать \[P(X_1 + X_2 > 0.8 \mid X_3 < 0.1)\]

n_obs <- 10^6
x_1 <- runif(n_obs, min = 0, max = 2)
x_2 <- runif(n_obs, min = 0, max = 2)
x_3 <- runif(n_obs, min = 0, max = 2)

uslovie <- x_3 < 0.1
x_1[1:5]
x_3[1:5]
uslovie[1:5]
x_1[uslovie][1:5]
## [1] 1.7260578 1.9716828 0.9364001 0.5078184 1.9227958
## [1] 0.7036184 1.2518040 0.5209894 1.4464742 1.6958722
## [1] FALSE FALSE FALSE FALSE FALSE
## [1] 0.2269900 1.5824290 0.5807984 0.9787026 0.6801348

Зачем всё это надо?

  • \(X_1, X_2, X_3 \sim U[0;2]\), независимые

  • Хотим знать \[P(X_1 + X_2 > 0.8 \mid X_3 < 0.1)\]

n_obs <- 10^6
x_1 <- runif(n_obs, min = 0, max = 2)
x_2 <- runif(n_obs, min = 0, max = 2)
x_3 <- runif(n_obs, min = 0, max = 2)

uslovie <- x_3 < 0.1

success <- x_1[uslovie] + x_2[uslovie] > 0.8

sum(success)/n_obs
## [1] 0.046189

Генерация выборок

С такой же лёгкостью можно генерировать любые выборки

sample(1:10, size = 8)  # выборка без повторений

sample(1:10, size = 8, replace = TRUE) # с повторениями
## [1]  9  2  1  5 10  8  7  3
## [1] 1 6 3 2 9 6 3 2
# неправильная монетка 
sample(c("Орёл", "Решка"), size = 5, replace = TRUE, prob = c(0.3, 0.7))
## [1] "Орёл"  "Решка" "Решка" "Решка" "Решка"

ЗБЧ и монетка

x <-  sample(c("Орёл", "Решка"), size = 10^6, replace = TRUE, prob = c(0.3, 0.7))

sum(x == 'Орёл')/length(x)
## [1] 0.300226

Центральная предельная теорема (ЦПТ)

ЦПТ

  • При определённых условиях сумма достаточно большого числа случайных величин имеет распределение близкое к нормальному

  • Главное: чтобы случайные величины были похожи и не было такого, что одна резко выделяется на фоне остальных

  • Есть много разных ЦПТ с разными условиями

Классическая формулировка ЦПТ

Пусть \(X_1, \ldots, X_n, \ldots\) — последовательность независимых одинаковых случайных велчин с конечным вторым моментом \(E(X_i^2) < \infty\). Тогда


\[ \frac{(X_1 + \ldots + X_n) - n \cdot E(X_1)}{\sqrt{n \cdot Var(X_1)}} \overset{d}{\to} N(0,1) \]

ЦПТ и равномерное распределение

  • Пусть \(X \sim U[-1;1]\), пусть \(Y = X_1 + \ldots + X_n\)

Доска Гальтона

ЦПТ на пальцах

  • \(Y\) — время прихода на первую пару

ЦПТ на пальцах

  • \(Y\) — время прихода Саши на первую пару

  • На Сашу прыгнул кот и он проснулся пораньше, ускорение на \(X_1\)

  • Пока готовил завтрак, убежало молоко, задержка на \(X_2\)

  • Быстро приехал автобус, ускорение на \(X_3\)

  • Встал в неожиданную пробку, задержка на \(X_4\)

ЦПТ на пальцах

  • Ни одна из величин не выделяется!

ЦПТ и ЗБЧ как скорости

В ЗБЧ:

\[ \frac{X_1 + \ldots + X_n - n E(X_1)}{n} \overset{p}{\to} 0 \]


В ЦПТ та же дробь домножается на \(\sqrt{n}\), это замедляет сходимость и мы приходим к более интересному результату:

\[ \sqrt{n} \cdot \frac{X_1 + \ldots + X_n - n E(X_1)}{n} \overset{d}{\to} N(0, Var(X_1)) \]

Крайнестан и среднестан

  • ЦПТ и ЗБЧ работают в среднестане

  • А что, если какая-то одна случайная величина выбивается?

  • Тогда мы перемещаемся из среднестана в крайнестан и сталкиваемся с проблемой тяжёлых хвостов

  • О тяжёлых хвостах мы ещё поговорим, они часто выскакивают в мире финансов

Фальсификации

Примеры случайностей из среднестана

  • Хочу примеры случайных величин, которые складываются из суммы других и из-за этого имеют нормальное распределение!

Примеры случайностей из среднестана

  • Хочу примеры случайных величин, которые складываются из суммы других и из-за этого имеют нормальное распределение!

  • Рост человека, вес человека

  • Погрешность весов

  • Выборочное среднее, посчитанное по выборке объёма \(100\)

  • Умственные способности человека, его балл человека за ЕГЭ

Результаты ЕГЭ в Польше

Результаты ЕГЭ в Польше

  • Подозрительный пик в районе проходного балла

  • Подозрительный пик на \(100\) баллах

Явка на выборы

  • Явка на выборы вроде бы находится в среднестане, выборы довольно масштабное явление

Выборы в европейских странах

Последние выборы в Германии и Франции

Выбры в ещё парочке стран

Выборы в России

  • Странное распределение явки, пики в круглых числах: \(55\), \(60\), \(65\), \(70\)

Выборы в Москве

В Москве с распределением явки всё окей

Выборы в Кемерово

В Кемерово наблюдается пила Чурова

??????????

  • Голосование в России противоречит ЦПТ и результатам электоральных исследований по другим странам

  • На волеизъявление людей что-то повлияло

Ответ ЦИК

  • ЦИК при Чурове в ответ на критику написал свою научную работу

  • Множество избирателей неоднородно, в сёлах массово сильны провластные взгляды, нужно детализировать подмножества избирателей

Зависимость от явки

  • Главный вопрос в зависимости результатов одной партии от явки

Корректировка результатов голосования

Выводы

  • Не факт, что явка на выборы должна быть куполообразной в неоднородных странах

  • Чисто теоретически доля кандидата может зависеть от явки, но это неточно

  • Тем не менее, декомпозиция по кандидатам и участкам с автоподсчётом голосов могут являться косвенным доказательством наличия вбросов

Полезные ссылки

Случайности бывают разными

Распределение Бернулли

Биномиальное распределение

Распределение Пуассона

Любое дискретное распределение

Вырожденное распределение

Экспоненциальное распределение

Равномерное распределение

Связь распределений

Безумные рекомендации